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We examine the dynamical behavior of recently introduced bubbles in asymptotically flat, five- 
dimensional spacetimes. Using numerical methods, we find that even bubbles that initially start 
expanding eventually collapse to a Schwarzschild-Tangherlini black hole. 
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I. INTRODUCTION 

o 

"Bubbles of nothing" in higher dimensional spacetimes have been the subject of significant attention in recent years. 
£SJ . They have played a key role in understanding the phase space of black hole spacetimes in Kaluza-Klein scenarios [3, 0] , 
have surfaced as a possible way around the black hole information paradox [3| and as mediators of non-perturbative 
instabilities in AdS/CFT contexts [J], been discussed in connections with orbifold decays in AdS spacetimes Q and 
have been the subject of studies of bubble-bubble collisions 0, 0)- However, essentially all higher dimensional bubble 
studies have been presented within Kaluza-Klein scenarios as no known non-local bubble solutions outside them were 
. available. Recent work by Copsey 0, 0] presents the first examples of bubbles in asymptotically flat and AdS cases. 
In particular, Ref. [9( provides data, at a moment of time symmetry, where surfaces of constant radius r are squashed 
r 3 1 , three-spheres with the circumference of one of the circles converging to zero when r approaches some fixed positive 
value while these surfaces are metric three-spheres in the asymptotic region r — * oo. By judiciously choosing the 
freely available geometrical variables, a "bubble of nothing" can be defined. Since Copsey's construction gives data at 
• a moment of time a complete picture of the spacetime is not available. While several important observations can be 
\ drawn at the hypersurface where this data is known, like the initial growth-rate of the bubble, key issues can only be 
l— ~~ '■ addressed by examining the full evolution of the spacetime subject to the provided initial data. Among the questions 
one would like to answer are 

• What is the bubble's behavior in time? In particular, if a bubble begins expanding, does it keep on expanding? 
*0 ! If the expansion were exponentially fast as in the case of Kaluza-Klein bubbles [101 ] , the spacetime might have 

the ingredients for realizing a possible violation of cosmic censorship. This would occur if the bubble expanded 
, so as to meet (a portion of) future null infinity at a finite affinc time of the generators of I + . 

\Q • • If the bubble were to reverse its expansion rate, what is the end behavior? Does the bubble collapse? If so, does 

f^*) ' it form a black hole? 

; 

• If a black hole is formed, is the bubble size when the horizon forms already at the string/Planck scales? In 
this case the classical evolution certainly could not be trusted and string effects should be taken into account 

. ' to reveal the bubble's final fate. 
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In this work we examine the evolution of some of the bubble initial data presented by Copsey and answer the above 
questions. 

Our article is organized as follows. In section [II] we review Copsey's initial data representing exact bubble solutions 
at a moment of time symmetry. The future of this data is then obtained by solving the Cauchy problem by numerical 
methods. The system of evolution and constraint equations, the coordinate choices and boundary conditions used 
to evolve this initial data set as well as the numerical implementation are discussed in section IIIII In section IIVI we 
present the results of the evolutions carried out and we end in section |V] with a brief discussion and conclusions drawn 
from this work. Technicalities like the computation of the curvature tensor and curvature invariants and a summary 
of the Schwarzschild-Tangherlini solution can be found in the appendices. 



II. PROTOTYPE METRIC AND INITIAL DATA 



Recently, Copsey [9] introduced bubbles which are outside the traditionally Kaluza-Klein framework. These bubbles, 
for the asymptotically flat case, and at a moment of time-symmetry, are described by a t — const slice of a metric of 
the form 



ds z 



-dt z 



dr 2 
W{r) 



r 

T 



H(r) {dz + cos dip) 2 + dtf 2 + sin 2 <Q dip 1 



(1) 
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where W and H are smooth functions of r which converge to 1 as r — > oo and where z) € [0, ir) x [0, 27r) x [0, 4ir) 

are Euler angles on the three-sphere S 3 . If W = H = 1, this metric describes five-dimensional Minkowski space. This 
can be seen by introducing the coordinates 



x + iy = r e 2 ^ z ^ sin 



X + iY = re^ {z+v) cos 



in terms of which the metric ([T|) with W = H = 1 assumes the form 

ds 2 = -dt 2 + d.T 2 + rfy 2 + dX 2 + dY 2 . 

More generally, using 

dr = riidx 1 , dz + cos d dip — rriidx 1 , 
where (x l ) = (x, y, X, Y), (m) — (x,y, X. Y) jr and (m,) = {—y, x, — Y, X)/r, we find that 

1 - W 



ds 2 = —dt 2 + h ti dx l dx j = —dt 2 



W 



(H — Vjmirrij 



dx l dx J . 



Therefore, the metric ([T]) is asymptotically flat if the functions W — 1 and H — 1 decay to zero fast enough as r — > oo. 
In particular, if W = 1 - Ci/r 2 + 0(r~ 3 ), H = 1- C 2 /r 2 + 0(r~ 3 ), iT = 2C 2 /r 3 + 0(r" 4 ) we obtain a finite ADM 
mass 
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1 f 4 

— lim / V (ft^ - djhii) dSj = I (3d - C 2 ), 

lD7T ' — >oo J • 4 7^ 8 



where S 13 denotes the three-sphere with radius r and rfS^ = njdS is the area element on S 3 . 

Copsey's construction also assumes that the functions W and H both have a single root at some ro > and are 
both strictly positive for r > tq. This means that the circumference of the circles determined by the orbits of the 
Killing field d z shrinks to zero at r = tq. In order to understand the geometry near r = ro, we replace r by the new 
radial coordinate 



R = dr 2 -r 2 , 



r > r . 



Then, the metric ([T]) can be rewritten as 

R 2 dR 2 



ds z = -dt 2 



R 2 



(R 2 +r 2 )W 4 
Because of the assumptions on W and H, they have the form 



H (dz + cos & dip) + dd 2 + sin 2 d dip' 



(2) 



W 



H 



2a z 



2(3' 



ro 



ro 
r - r 



O 



O 



ro 



ro 



ro 



r-r 
ro 



= fi 2 







o 



near R = 0, where a and (3 are two strictly positive constants. Therefore, as R tends to zero, the metric has the form 



dR 2 

ds 2 ~ -dt 2 + — + tL-R 2 (dz + costfdip) 2 + '-^ [dd 2 + sin 2 d dip 2 ] . 

Since z has period 47r, there is a conical singularity at R — unless a[3 = 1. If a(3 = 1, we may replace R and z with 
the Cartesian coordinates (u,v), which are defined by 

R 



a 2 



u + iv 

V 

In terms of these, the metric ([2]) assumes the form 

r 2 , + a(u 2 + v 2 ) 



iz/2 



ds 2 



-dt 2 + du 2 + dv 2 + 



[d$ 2 + sin 2 dip 2 



+ f(u 2 + v 2 )(udu + vdv) 2 +h(u 2 +v 2 )[2(udv-vdu) + (u 2 +v 2 )cos'& dtp] , 
where / and h are smooth functions of u 2 + v 2 . Therefore, if a/3 = 1, we obtain a smooth, regular asymptotically flat 



manifold with topology I 
Killing field d z vanishes. 



x S . At each fixed time t, the bubble is determined by the two-sphere at which the 
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A. Initial data 



Data at the moment of time symmetry automatically solves the momentum constraint and one is left with having 
to satisfy the Hamiltonian constraint only. For a ( = const section of a metric of the form ([T]), the Hamiltonian 
constraint in vacuum yields 



IT 



4 if' 
r~H 



2 H 2 



6 





lif' 3" 


w + 






Ye + r 



w 



(if - 4) = Q, 



(3) 



where a prime denotes differentiation with respect to r. A convenient way to solve this equation is to freely specify 
H and integrate the resulting linear equation for W . In order to analyze this, let H be an arbitrary smooth function 
with the following properties: H has a single root at some point r = r > 0, H(r) > for all r > r and H(r) — 
1 - C 2 /r 2 + 0(r~ 3 ), H'(r) = 2C 2 /r 3 + <3(r~ 4 ) and H"(r) = -6C 2 /r 4 + <3(r~ 5 ) for r — s- oo. Furthermore, we assume 
that the expansion along outgoing null radial geodesies, which is proportional to 6H + rH', is strictly positive for all 
r > ro. This condition ensures the initial data does not contain any apparent horizons Q. 

For the following, it is convenient to replace the function W by a new function £, defined by W(r) = H(r)(\ + ((r)), 
and to introduce the dimensionless compactified coordinate s := ro/r which varies from to 1. Then, Eq. Q reads 

^(6if - sH s )( s = (s 2 H ss - 5sH s + 6H)C + [s 2 H ss - 5sH s + 8(H - 1)] , 



where the subscript s denotes differentiation with respect to s. The general solution to this equation is 

C(*)=Co(«) [C + J(s)], 



(4) 



where Cq is a constant and 



CoOO 



(6H - sH s ) 2 



s 

J(s) = 2 J{<oH- tH t 



t 2 H tt -5tH t + 8(H-1) e -2i(r) dT 



with 



/(*) 



H T dr 



6H - tH t 



Because of the above assumptions on H, these integrals are well-defined and converge to a finite value as s — > 1. In 
particular, notice that the asymptotic behavior on H implies that t 2 H tt — 5tH t + 8(H — 1) = 0(t 3 ) so that the 
integrand in the expression for J(s) is well-defined at t = 0. Furthermore, we see that CO) = 0(s 2 ) for small s > 0, 
hence the function W has the same asymptotic behavior as H and possesses a root at r = ro- The constant Co is 
determined from the requirement of the absence of a conical singularity at r = ro, H S (1)W S (1) — 4, which yields 

C =[4-tf s (l) 2 ] e- 2 '«-J(l). 

Finally, we notice that without further restriction on H there is no guarantee that W is positive for all r > r$. 
However, the condition H < 4 is sufficient to guarantee positivity of W for r > tq. Indeed, Eq. ([3]) and the condition 
6H + rH' > then imply that W'(r±) > if r\ is a zero of W which shows that the sign of W cannot change from 
positive to negative when r increases. 

In order to understand the initial dynamics of the bubble, we compute the initial acceleration of the bubble area A 
with respect to proper time r at the bubble. In terms of the general metric (IA1[) the bubble's area and proper time 
are given by 



A = Aire 



2c I 



eh = e d dt\ 



respectively. Taking the second derivative with respect to proper time, using the evolution equation R33 
appendix [5] and the initial values 



in 



R = r, 



■logOn 



log 



log(if), 



log 



4 



in (|Aip , we obtain 



d 2 A 

dT 2 



-2irW s {l) - 8tt 



-8tt 



H S {1) 



(5) 



Therefore, the bubble starts expanding if and only if — 1 < H s (l) < 0. Finally, the ADM mass is 



M A 



7rrn 



DM 



Co 



Co 
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Some examples of solutions, and the ones we will employ in our simulations are 
1. The choice W — H simplifies Eq. © to the following linear equation for H, 



which has the general solution 



H" 



H(r) 
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-H' + ^(H-1) = 0, 



1 - a ( — 



(6) 



with a,Q + ai — 1. Notice that H(ro) = and that H can also be written in the form 



H(r) = 



r 



1 



ai — 

r 



from which we see that a necessary and sufficient condition for H to have a single root at ro and to be strictly 
positive for r > ro is that ai > —1, or, equivalently, ao < 2. Furthermore, H converges to one as r — > oo. In order 
to verify the regularity condition at r = ro, we compute 2a 2 = roH'(ro) — 2(1 + ai). Hence, a 2 = 1 = 1 + ai 
and the regularity condition a/3 = 1 forces the choice ao = 1, a\ = 0. 

2. More general initial data can be obtained by employing the function H given in Eq. © but relaxing the 
condition W = H. In this case, the integrals in Eq. (|U) can be performed analytically [9|] with the result 



ai(^) 2 + a + 6c 



with 



6n = 



3ai 



ai + a + foo 

ai(ai + 2) 
(1 + ai) 2 ' 



(7) 



If ai = 0, £ = and we recover the previous solution with W — H. However, the more general family of 
solutions given by H as in Eq. ((6J and W = H(\ + £) with ( given in Eqs. ((7]) allows for initially collapsing 
and expanding bubbles [9] since the initial acceleration is 



d 2 A „ 2a - 3 
——r- = 4tt 

dr 2 2 — ao 



(8) 



Therefore, the bubble initially expands for 3/2 < ao < 2 and initially collapses for ao < 3/2. The parameter ao 
determines the ADM mass through 



M A 



DM 



2a + 3 



ai(ai + 2) 
(1 + ai) 2 



a - b 



ai + a - b 



ao + b 



ai + ao + bo 



2b 



3. A simple and interesting example is given by the choice H(s) = (1 — s 2 ) 2 + es 2 (l — s), where < e < 1. 
This choice satisfies all the required assumptions on H, and < H(s) < 1 + s for < s < 1. Furthermore, 
H s (l) — —e, so one can construct bubbles with arbitrarily large initial acceleration. 

4. The following ansatz is also suggested in [9(, 



H(s) = (l-.s 2 ) 4 + 4es 2 (l-s) 



4s 2 (l-s) 2 c 2 



(9) 



which leads to initially hyperexpanding bubbles. 
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III. EVOLUTION 



Here, we discuss our method for obtaining the time evolution of the time-symmetric bubble configurations discussed 
in the previous section. Motivated by the form ^ of the prototype metric, we perform the following rescaling in the 
general line element (|A1[) 



d. 



I>+ !<-:> ( | 



log 



(10) 



The metric now takes the form 



ds 2 



e 2d dt 2 + e 2a dR 



2a, 



e 2b R 2 {dz + costf d<p) z + e 2c (R 2 + 4)dQ 



R>0. 



(11) 



The functions a, b, c and d which depend only on t and R must satisfy the following boundary conditions. As R — ► oo 
we require that these functions converge to zero fast enough for quantities like the ADM mass to be defined. Near 
R = we require that these functions are smooth and satisfy the conditions a'(t, 0) = b'(t, 0) = c'(i, 0) = d'(t, 0) = 0. 
Furthermore, we impose the condition a(t, 0) — b(t, 0) = which ensures that there is no conical singularity at R = 0. 



A. Gauge conditions 



We find it convenient to impose the following family of gauge conditions on the logarithm of the lapse d, 

d = a + \{b + 2c), (12) 

where A is a fixed parameter. For A = 0, this condition implies that the two-metric —e dt 2 + e 2a dR 2 is in the 
conformal flat gauge. As we will see, the principal part of the evolution equations is governed by the d'Alembertian 
with respect to this metric. Since the two-dimensional d'Alembertian operator is conformally covariant, the resulting 
evolution equations are semi-linear. In particular, this implies that the characteristic speeds do not depend on the 
solution that is being evolved, so there cannot be shock formation due to the crossing of characteristics. When 
A = 1, the gauge condition (|12p is strongly related to the densitized lapse condition often encountered in hyperbolic 
formulations of Einstein's field equations. Indeed, the square root of the determinant of the four metric belonging to 
(fTTjl is e a+b+2c R(R 2 + 7"o) sini9/8, so (|12[) sets the lapse e d equal to the square root of the determinant of the four 
metric divided by the factor R(R 2 + Tq) sin$/8 which is singular at the bubble, at the poles d = 0, ir and at R = oo. 
Since this gauge condition is essentially the time harmonic gauge condition integrated in time, it is convenient for its 
singularity avoidance behavior. In our simulations below, we use both choices A = and A = 1. 



B. Evolution equations 

In the gauge (fT2")) . the evolution equations can be written as a coupled system of three wave equations for the fields 
a, b and c which are obtained from the Einstein vacuum equations R%\ — (A + l)Goo = 0, R22 = and i? 33 = 0, 
respectively. The resulting system has the form 



= 2A(b+2c) 



-Mu' 
R 



+ F(R, u, it, u'), 



where u := (a, b, c) T , M is the constant 3x3 matrix 



A -2A -2(A + 1) 
\l =■ I A + 2 2(A + 1) 
1 

and J 7 is a nonlinear function of R, u, u and v! which is regular at R = 0. Provided that u is smooth enough, the 
evolution equations are regular at R = since then the boundary condition u'(t,0) = implies that 



lim — Mu' = Mu" 
R^O R 



R=0 
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For the following, we find it convenient to replace a, b and c by the linear combinations A = a + Xb + 2(A + l)c, 
B = b + 2c and C = c which diagonalize M. Using the expressions in the appendix and taking into account the 
rescaling (fTU| we obtain 



.4 



D 



1 



A„2AB^n' 



2XB 



-3(A- 1)C" 2 + 



2i? 



i? 2 



(AA' + 2XB' - 3(A - 1)C") + 



2(l + A)r 2 + (3A + l)i? 2 



+ 2A.4.B - A(A + l)B 2 - 3(A + 1)C 2 + (A + 3)Vb - (A + l)Vi , 



i? A + 2 



a (A-l)B 



R X+2 e (X+l)B B , 



J2XB 



R 2 



[(A + 2)RB' + 3] + (A - \)B 2 + 2V - 2Vi , 



R V 



„2AB 



r 2 + R 2 



[(A + l)RB' + 2RC' + 2] + (A - l)BC + 2V - V x , 



(13) 
(14) 

(15) 



where V Q = R 2 (r 2 + R 2 )~~ 2 e 2 ( A+B ~ 6C *> and V\ = 4(r^ + R 2 )~ 1 e 2t - A ~ 3C '\ These equations are regular at R = provided 
the boundary conditions A' = B' = C = at R = are satisfied and the fields are smooth enough. Eqs. (|13I14I15[) 
are subject to the Hamiltonian constraint Ti. = and to the momentum constraint M. = 0, where 



H = -e 2d G 00 

P (A-1)B r 



R X+2 e (X+l)B B , 



2XB 



R 



A' + A'B' - 3C*' : 



2R 



R 2 



(A' - (A + 1)B' - 3C") 



4r 2 + 3i? 2 



(r 2 +i? 2 ) 2 



AB + AS 2 + 3G 2 +V -Vi, 



M 



e 2d R i 



„AB 



B' - AS' - B {A 1 - (A + 1)5') + 6GC" - i (A - (A + l)s) - 2 2J? (A - XB - 3(7 



Notice that 



lim RM = - e XB (A-(X+ 1)b\ 



R=0 



(16) 
(17) 

(18) 



hence the satisfaction of the momentum constraint implies the condition a(t, 0) — b(t, 0) = A(t, 0) — (A + l)B(t, 0) = 
for the avoidance of the conical singularity at R = 0. In the next subsection, we show that the evolution equations 
and the regularity conditions A'(t,0) = B'(t,0) = C'(t, 0) = imply that the constraints Ti = A4 = are satisfied 
everywhere and at each time if satisfied initially. 



C. Propagation of the constraints 



As a consequence of the twice contracted Bianchi identities (sec appendix [XJ and the evolution equations which 
imply Gn — Goo, G22 = G33 = G44 = — AGoo we find that the constraint variables Ti and M. obey the following 
linear evolution system 



U 



(3A - l)BH + e 



XB 



M'+\B' 



r 2 + 3R 2 
R{r 2 + R 2 ) 



M 



M = (2X-l)BM 



.AS 



n 1 + (1 + a) b' 



r 2 + 3R 2 
R(r 2 + R 2 ) 



This system can be simplified by introducing the rescaled variables 

H := R{rl + R 2 )e B Ti, M :— i?(r 2 + R 2 )e B M. 
In terms of these, the constraint propagation system reads 

Ti = 3XBH + e XB M', 



M = 2XBM 



XB 



H' + X[B' 



r 2 + 3R 2 
R(r 2 + R 2 ) 



n 



n 



(19) 
(20) 



(21) 
(22) 
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Notice that the regularity conditions A'(i,0) = B'(t,0) = C"(t,0) = at R = imply that |m « = and that 
lim M. exists. Therefore, Eqs. (|21l22p are regular at R — 0. Defining the following "energy" norm for the constraint 

R — >0 

variables 

oo 

£{t) := i J (ft 2 + M 2 ) R x (r 2 + R 2 ) x dR, 
o 

taking a time derivative and using Eqs. (|21I22[) we obtain 



dt 



3C 



^-£ = R x {r 2 + R 2 ) x e XB HM + A / B ( 3H Z + 2M 2 ) i? A (ro + R' 2 ) A dR. (23) 



i?=0 



n 



The boundary term on the right-hand side vanishes because of the regularity conditions at R = and under the 
assumption that all fields fall off sufficiently fast as R — > oo. If the quantity B is bounded, we can estimate the 
integral term on the right-hand side by a constant K times the energy norm £ , and we obtain an estimate of the form 

S(t) < e^ Kt £(0). 

This shows that if the constraints are satisfied initially, they are also satisfied for alH > for which a smooth enough 
solution to the evolution equations exists. More generally, this inequality shows that convergence of £(0) to zero 
implies that £ (t) also converges to zero for each such t. This is important for the numerical implementation below 
where the initial data does not satisfy the constraints exactly due to truncation errors, but where the initial constraint 
violation converges to zero with increasing resolution. Finally, we notice that in the particular gauge where A = the 
"energy" norm £ cannot grow in time. 



D. Outer boundary conditions 



In our numerical implementation below, we truncate the domain at some large radius R m ax 2> ?~o and impose the 
following boundary conditions at R = R ma x- First, we enforce the momentum constraint M\ R=R — 0, which 
yields the boundary condition 



B' 



R=R„ 



AB 1 + B(A' - (A + l)B') - 6CC + i (A - (A + l)B) + ^ ™ R2 (A-XB- 3c) 



(24) 



R=R„ 



When estimating the constraint errors, oo has to be replaced by R ma x in the expression for £(t) above. One then 
obtains the same equality as in Eq. (|23|) but with oo replaced by R m ax on the right-hand side of that equation. The 
enforcing of the momentum constraint then guarantees that the boundary term vanishes. Therefore, we conclude as 
before that the constraint errors converge to zero at each fixed t > if the initial constraint errors converge to zero. 
Finally, we impose a Sommerfeld-like boundary condition on the fields A and C, 



A + e A' 



R=R n 



0, 



C 



c" 



R=R„ 



0. 



(25) 



These two boundary conditions set the incoming characteristic fields corresponding to the evolution equations for 
A and C to zero, and thus do absorb high-frequency outgoing waves. However, these boundary conditions do yield 
reflections for low-frequency waves propagating towards large R. In our simulations below, we choose Rmax to be 
large enough so that such reflections do not influence the region where physics is extracted. 



E. Numerical implementation 

In order to numerically implement the evolution system (I13I14I15P we find it convenient to rewrite it as a first order 
symmetric hyperbolic system by introducing the variables tta ■— A, ttb ■= B, ttc '■= C and ipA ■— A', ips '■= B', 
ipe ■= C. The resulting system has the form 

U = 7T, (26) 

i = e MB r 1 A W ) + f( M ,^;4 (27) 

* = ^ (28) 
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where u = (A,B,C), ir = (tt a , n B , n c ), ip = {ip A ,ipB,ipc), E = dmg (R x e 2XB , R x+2 e^+^ B , Re^ +1 ^ B ) and 
F(u, it, ip; R) is a nonlinear term that can be read off from Eqs. (|13I14I15[) . This system is then discretized by 
the method of lines. Let us start with the description of our spatial discretization. 

We consider first a uniform grid Rj = jAR, j — 0, 1, 2, ...N, where AR > is a fixed mesh size, and replace the 
functions u, it and i/j with gridfunctions Uj, ttj and tpj, j = 0, 1, 2, ...N, respectively. The evolution system (|26I27I28[) 
is then approximated by the semi-discrete system 



3 > 

e 



(29) 



^'E^D^.+Fiu^^R,), j = l,2,..JV, (30) 



u 3 — y^n/j , 

where D denotes the finite differencing operator 



i>j = (Dir) 3 , (31) 



^j-i „• -I o AT — 1 

2AR ' J ~ < z > 1 ' 
AR ' J iv ' 



which is second order accurate in the interior points j = 1,2,..., TV — 1 and first order accurate at the exterior point 
j = N. At the inner point j = 0, we use the regularity condition ip(t, 0) = and assume smoothness of ip, which 
imply 



n-l ^ tT^ 1\ A dip 



lim E- 1 — (Eip) =A-j§ , A = diag(A + l,A + 3,2) 
-R^o dR dR 



by de L'Hopital's rule. We then evolve mq, ttq and tpo according to 



UQ = 7T , (32) 

= e 2XB °A^+F(u, tt.^O), (33) 

= 0. (34) 
In order to impose the outer boundary conditions (|24l25p we replace (i}>b)n by 



BIN 



1 27? 

7TA-0B + ttb(^a - (A + 1)V>b) - O^c^c + -5 (vr^ — (A + 1)tt b ) + 2 (tta - A7r B - 37r c ) 

-ft + -ft 



JV 



and apply Olsson's projection method [ll[ in order to incorporate the Sommerfeld-type conditions (|25|) . This method 
consists in applying the projector 



as 

2 \ XB 1 



to the two-vectors (tta,iI>a)n and (ttc, iPc)n- The projected fields then satisfy 

(tt a ) w + e AB " (ip A )N = 0, (ttc)jv + e AB " {^ c )n = 0, 



which agrees with (|25|) in the high-frequency approximation. 

In our simulations, we choose units in which r*o = 1. In order to improve the accuracy of our code, we replace the 
radial coordinate R by a computational coordinate x with 

x 

R = 1 j x ) < a; < a; ma:c , (35) 

where < k < 1 is a "stretching" parameter. A uniform grid Xj — jAx, j — 0, 1,2, ...N in x is then used, where 
Ax = x rnax /N. Since dR = dx/(l ~ kx / x max ) , a uniform grid in x implies (for k > 0) a non-uniform grid for R, 
with a smaller mesh size near the bubble than near the outer boundary. Since the fields have their largest gradient 
near the bubble, this non-uniformity in R helps better resolve the dynamics. Typically, we use the values x max = 10 
and k — 0.8 in our simulations which implies that R m ax is fifty times as large as tq. 

Next, we use a third order Runge-Kutta algorithm in order to perform the time integration. In our simulations, we 
use 0.2 as a value for the Courant-Friedrichs-Levy factor At/ Ax. In [l2| we have shown for a related model problem 
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consisting of the reduced wave equation for spherically symmetric solutions that this discretization leads to a stable 
and convergent scheme. 

Finally, initial data for u, ir and ip is obtained by comparing the metrics ((21) and (fTTjl . This yields 
A = a + Xb, B = b, 
and 7Ta = ttb = 7rc = 0, where 

a = -Uog{W/R 2 ) 

and where we recall that W/R 2 and H/R 2 converge to a finite value as R — > 0. 

IV. RESULTS 

The numerical evolution of the initial data discussed in section |TT] suggest that the bubble eventually collapses to 
the Schwarzschild-Tangherlini solution, regardless whether or not the bubble is initially expanding or collapsing. This 
is shown in detail in the following subsections. 

A. Evolution of the bubble area 

In Fig. [1] we show the typical dynamical behavior of the bubble for type 2 and 4 initial data. More specifically, we 
show the bubble's area as a function of proper time at the bubble. It can be seen that even for type 4 initial data 
configurations with a violent initial expansion the bubble finally collapses in a finite proper time. We also explored 
a wider range of values for the parameters ao and e, also evolved type 3 initial data and found that the result is the 
same in all cases: the bubble area shrinks to zero after a finite proper time. 

The numerical parameters used in these calculations are: R ma x = 50, k = 0.8, N = 1000. In the case of the type 
2 data, we use the gauge condition ([12"]) with A = 1. Because of the singularity avoidance property of this gauge 
condition, we are able to continue the evolution for some time after the bubble area has shrinked to zero and are 
able to follow the evolution for long enough time to observe the formation of an apparent horizon and its settling 
down to an equilibrium black hole (see next subsection). However, for the evolution of type 4 data we observe the 
development of large gradients in the derivatives of the fields, similar to those found in a shock formation situation. 
This happens before the turning point in the evolution of the bubble area and we are unable to continue the evolution. 
For this reason we choose the gauge condition (fT2f with A = instead. As mentioned above the evolution equations 
are semilinear in this case, and thus there cannot be shocks due to crossing of characteristics. 

B. Formation of an apparent horizon 

In order to provide evidence for the bubble settling down to the Schwarzschild-Tangherlini solution, we first check 
whether an apparent horizon forms. This is done by calculating the expansion of the outgoing null normal defined as 
l a = n a + s a with n a the future-directed timelike normal to the t = const hypersurface and s a the outward pointing 
unit normal to the surfaces R = const which is orthogonal to n a . For the t = const surfaces of the metric (|A1[) the 
expansion reduces to 

9 = e- a (b' + 2c') + e- d (b + 2c). 

In terms of our main evolution variables (i.e. after the rescaling of the lower case variables, and then the redefinition 
in terms of the upper case ones), this can be reexpressed as 

e A-2C 6 = e XB ( B , 3R2 + r \ ) + B . 

\ R(R2 + r )J 



C = 0, 



da 
dR 



db 
'dR 



1pB 



db_ 

dR ' 



= 0, 



1 



■log(r 2 +R 2 



\og(H/R 2 ) + -log(r 2 



The apparent horizon is determined by the outermost surface defined by 9 — 0. Since the factor e A 2C is positive, it 
can be ignored for our purposes. 
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FIG. 1: We show the bubble area (normalized by its initial area) versus proper time evaluated at the bubble. Left: evolution 
of initial data of type 2 with several choices for the initial parameter ao. Right: evolution of type 4 initial data for C2 = 1 and 
different values of the parameter e. The inset shows an initial violent expansion of the bubble area. Nevertheless, this behavior 
changes and the bubble starts decelerating. In all cases the area shrinks to zero after a finite proper time. 

After detecting an apparent horizon we evaluate the curvature invariants In , Jn and Kn defined in appendix [B] 
where we replace Teh by the apparent horizon radius tah defined as the areal radius of the t — const cross sections 
of the apparent horizon, 

r AH :=[R(r%+R*)e B } 1/3 . (36) 
In order to check that the horizon is 50(4) symmetric, we also compute the quantity 

_ e 26 _ e 2c _ R 2 e 2(B-2C) _ {R 2 + ^yC 
M ''~ 6 2b + e 2c ~ R2 e 2(B-2C) + (#2 + r 2) e 2C ' 

on the apparent horizon. This quantity measures how close a given three-surface t — const, R — const, is of being a 
metric S 3 . In particular, a = if this surface is a metric (unsquashed) S* 3 . 

The time evolution of the invariants In, Jn, Kn and /x are shown in Fig. [2] for two different values of the parameter 
ao for type 2 initial data. It can be seen that In, Jn and Kn converge to 1 and that \i converges to zero, which 
is a strong indication that the apparent horizon settles down to the event horizon of the Schwarzschild-Tangherlini 
solution. The numerical parameters used in these calculations are: R max = 50, k = 0.8, N = 4000, A = 1. For type 
4 initial data, where we need to use the gauge choice (fl~2"|) with A = in order to get to the collapse, we are able to 
detect the formation of an apparent horizon, but after that we soon have to stop the simulation because the resulting 
gauge condition is not singularity avoiding. In order to track the late time behavior of the invariants we switch from 
A = to A = 1 shortly after the bubble starts collapsing. We then observe the same qualitative features as for type 2 
initial data: an apparent horizon forms and the curvature invariants In, Jn and Kn converge to one and /z to zero. 
Therefore, it seems that in this case too the solution settles down to a Schwarzschild-Tangherlini black hole. 

As indicated, all cases explored where the bubbles begin expanding turn around and eventually collapse. This 
collapse is hidden however from external observers by the appearance of an event horizon before the bubble shrinks 
to zero size. An important issue to examine is the size of the bubble when a horizon forms as if it were significantly 
smaller than the bubble's initial size string effects would become relevant and the present classical analysis would not 
apply. To examine this behavior we consider the type 2 data family of initially expanding bubbles and obtain the size 
of the bubble when an apparent horizon first appears. The results (shown in table [J) indicate that the bubble only 
shrinks to about 1/3 of its initial size by the time the horizon forms. Therefore, a "classically-sized" bubble would 
seem to form a horizon while itself is still in the classical regime as suggested by its apparent horizon size. Certainly 
one should be careful with this argument since the apparent horizons are slicing-dependent. However, by monitoring 
the behavior of null rays one can also study the formation of an event horizon and the overall causal structure of 
the resulting spacctimc. To this end we trace past-directed ingoing null rays from a late-time slice of our simulation. 
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FIG. 2: The time evolution of the invariants In, Jn, Km and fx for type 2 initial data with two different values of ao. The 
apparent horizon forms at a time between t — 0.8 and 0.85 in the case with ao = 1.2 and at a time between t = 1.25 and 1.3 
for a = 1.99. 

This is shown in Fig. [3] for the type 2 initial data with ao = 1.99. As seen from the plot, the event horizon does not 
intersect the initial hypcrsurface but branches off the bubble shortly after the the bubble reaches its maximum size. 
We would like to mention that while the results discussed above, and the ones included in the table illustrate what is 
observed for the case of type 2 initial data, a similar behavior is obtained for the other data types discussed. 



ao 


1.55 


1.65 


1.75 


1.85 


1.95 


1.99 


Bubble size 


0.27 


0.27 


0.28 


0.29 


0.29 


0.29 



TABLE I: Size of the bubble when the apparent horizon first appears for different values of the initial data parameter ao- 




area! radius areal radius 



FIG. 3: Left: the past-directed ingoing null rays for the numerical spacetime obtained from type 2 initial data with ao = 1.99. 
The dotted line corresponds to the location of the apparent horizon and the dotted-dashed line indicates the location of the 
bubble area. Right: only a small subset of rays near the location of the apparent horizon is shown. This indicates that an 
event horizon forms at the bubble shortly after the bubble reaches its maximum size. As can be seen, the apparent horizon at 
late times t > 2.5 is a good approximation for the event horizon. Here, time t refers to coordinate time and not proper time 
as in Fig. [JJand the areal radius refers to the geometric radius of the two-spheres which is equal to e c \/ R 2 + 7q/2. Although 
not apparent from the plot, we have verified that the area of the S 3 cross sections of the event horizon does grow in time, as 
expected from the area theorem [l3j]. 

Further insight of the spacetime behavior can be obtained by comparing the ADM mass with the irreducible mass 
of the formed horizon. In Fig. [5] we compare the ADM mass with the apparent horizon mass of the final black 
hole which is defined as Mah = ^- r AH where Tah is the apparent horizon radius defined in Eq. (|5rj|) . The ADM 
mass is determined from the parameter an. in the initial data [9j. The difference between the two masses indicates 
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that significant radiation is produced during the collapse of the bubble. The numerical parameters used in these 
calculations are: R m ax = 50, k = 0.8, N — 4000, A = 1. 





FIG. 4: The Madm and horizon masses vs time are shown for the cases ao = 1.2 (left) and ao = 1.99 (right). 



C. Convergence checks 

In order to verify the convergence of our numerical approximation we performed several tests. First, we verify that 
the constraint errors become smaller as resolution is increased. This is shown in Figs. [S] and [5] for type 2 initial data 
with two different values of ag. In these figures, we compute the quantity CV = y + M-^ where 



Hoe '■= max 

l<j<N 



R i 



--H 



1 



3 > 



Mqo := max 

l<j<N 



Ri 



-.M 



1 



3 • 



As Figs. [5] and [5] indicate, the constraint errors show second order convergence to zero as resolution is increased. The 
numerical parameters used in these calculations are: R m ax — 50 and k = 0.8. 
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Case a = 1 .2 
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2 
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1.8 
1.7 
1.6 
1.5 



FIG. 5: Left: the quantity CV = y/H^ + M'^ vs time for various resolutions, for the value ao = 1.2. Right: the convergence 
factor CF C vs time. This factor is defined as CF C := log 2 (CV(iV = 8000)/CV{N = 16000)). 



Next, we also perform a self-convergence test for the field A. This is shown in Fig. [7] which indicates second order 
convergence. 

We also find convergence for type 3 and type 4 initial data although the corresponding convergence factors sometimes 
drop well below two in those cases, though remains above 1. This drop is due to the large gradients appearing in 
these solutions. 
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FIG. 6: The same as in the previous figure, but for the value do = 1.6. 



U 




"o = L2 ■ 

a n = 1.6 



FIG. 7: Convergence factor CF vs time for the case a = 1.6. This factor is defined as CF = log 2 (]|A(iV = 8000) - A(N = 
16000)||oo/||A(iV = 16000)- A(N = 32000)||oo), where || A||oo denotes the maximum value of A on the computational grid. The 
fact that CF is nearly equal two confirms second order convergence of our code. 



V. CONCLUSIONS 



In this work we examined the behavior of families of bubbles of nothing in fivc-dimcnsional asymptotically flat 
spacetimes. The simulations carried out reveal that independent of their initial behavior, these bubbles eventually 
collapse and form a black hole horizon around them. The formation of the horizon takes place at a time when the 
bubble's size is of the order of its initial size. Hence, starting with a 'classically-sized' bubble, the collapse and 
black hole formation occurs well above the string/Planck length. Thus, at least for the family considered, cosmic 
censorship is not violated and the overall process takes place at a classical level. The singularity that eventually forms 
in the spacetime is likely the same as that of the Tangherlini solution, thus its black hole singularity resolution at the 
quantum level should apply to this case as well. The behavior of initially expanding bubbles examined in this work 
is in sharp contrast with those with the same initial behavior in the Kaluza-Klein picture [lol ]. There the expansion 
continues and increases as time progresses and ultimately its expansion rate approaches light-speeds. This observation 
seems to indicate the presence of a fundamental underlying effect causing the turn-around in the asymptotically flat 
case. However, we want to stress that our analysis has not covered all asymptotically flat bubble solutions presented 
by Copsey and that we have not examined the AdS case. It would be extremely interesting to know whether the 
behavior found here is generic. 
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APPENDIX A: EXPRESSIONS FOR THE CURVATURE TENSOR 



In this appendix we compute the curvature and Ricci tensors belonging to the metric 



ds 2 = -e 2d dt 2 + e 2a dR 2 + e 2b {dz + costi dip) 2 + e 2c dn 



2b 



(Al) 



where d, a, b and c are functions of t and R only, and where dVt 2 denotes the standard metric on S 2 . In order to do 
so, it is convenient to introduce the following orthonormal basis of one-forms, 



6° = e d dt, e 1 = e a dR, 9 2 = e b (dz + cos dtp) , 6> 3 = e c di?, 
With respect to this, the curvature two-form Ojj = ^Rijkid A 6 l is 



6 i = e c sm$d<p. 



"01 — 


-k6° A0 1 , 






^02 = 


-0o 9" A e 2 


-Pi e 1 


A6 2 + (b - c)e b - 2c - d 9 3 A6» 4 , 


= 


-70 6° A 6 A 


-71 o 1 


a e A + 1(6 - c)e b - 2c - d e 2 A e A B0 B , 


n 12 = 


-/3i 0° a e 2 


-foe 1 


A 6 2 + (b' - c')e b - 2c - a 3 A 9 4 , 


QlA — 


-71 9° A e A 


-7 2 e 1 


A8 A + 1(6' - c')e b - 2c - a 9 2 A s AB B , 


&2A = 


-\(b-c)e b 


-2c-dg0 


A s A bO B - 1(6' - c!)e b - 2c - a 9 1 A e AB e B 


Q34 — 


(I - c)e b - 2c - 


- d e° a e 


2 + (b' - c')e b - 2c - a 9 1 A0 2 +S 2 6 3 A 9 4 . 



Here, A,B e {3, 4}, £34 = —£43 = 1, £33 = £44 = 0, a prime and a dot denote differentiation with respect to R and t, 
respectively, and 



K 

Po 

Pi 
02 
70 

7i 
72 

61 
62 



-2d 



-2d 



d + d(d-d) - e- 2a [d" + d'{d' - a')] , 



= e 



-a—d 



b + b(b- d) 

b' + b(b' -d') -ab' 

-2a \UI 1 1,1 (U 

-2d 



b'd', 



e ~[b" + b'(b' -a')] -e- 2d db, 



ii + c(c — d) 



— e 



2d; 

—2aJjf 



c'd', 



e - a - d [c' + c(c' -d')-a c'] , 
e~ 2a [c" + c'(c' -a')]-e- 2d ac, 



= e 



-2d 



1 



be- e- 2a b'c' + -e 
4 



26-4c 



e -2d t 2 _ £ -2a c ,2 _ - e 2b-4c + £ -2c 

4 



The corresponding components of the Ricci tensor are obtained from Rij = fl k i(ek, ej), where ej, j — 0, 1, 4, is the 
Fiinfbein corresponding to the orthonormal basis 9° , 1 , 6* 4 of one- forms. Explicitly, we find 



Rn 



-2d 



■d + d{a + b + 2c-d) - e- 2a [d" + b" + 2c" + d'(d' - a') + b'(b' - a') + 2c' (c' - a')} , (A2) 
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R 2 



-2d 



„-2d 



- 2a [b" + b'(b' + 2c' + d! - a')] + ^ e 



-a—d 



b + b(a + b + 2c - d) 
c+ c(a + b + 2c - d) 

b' + 2c' - d'(b + 2c) + b'{b -a) + 2c'{c - a) 



2fc-4c 



e- 2a \c" ' JIV 



R33 — R44 — e 

-R01 = Rio = 

Roo = 2Goo — R11 — R22 — R33 — R44 ■ 
Here, Gqo refers to the 00 components of the Einstein tensor, which is 



c'{b' + 2c' + d'- a')} - - e 2b - ic + e~ 2c , 



n — v- 2d 
Lroo — e 



a(b + 2c) + c(c + 2b) 



[b" + 2c" + {b' - a'){b' + 2c') + 3c' 2 ] - - c 2h ^ c + 



(A3) 

(A4) 

(A5) 
(A6) 

(A7) 



The remaining components of Rij are identically zero. 

Using the fact that the metric is 50(3) symmetric and homogeneous in the z direction, several curvature invariants 
can be defined. For this, consider the five- metric ds 2 in Eq. (IA1[) . the orbit three- metric 

ds 2 = -e 2d dt 2 + e 2a dR 2 + e 2b dz 2 

and the projection thereof on the spaces orthogonal to the Killing field d z 

ds 2 = -e 2d dt 2 + c 2a dR 2 . 

The Kretschmann invariants I = RijkiR l -' kl with respect to these metrics are 



I 



I = 4 [k 2 + p 2 - 2{i 2 + p 2 ] 



(A8) 
(A9) 



1 = 4 



k 2 + p 2 - 201 + p 2 + 2 7o 2 - 4 7l 2 + 2 7 | + 2<5 2 + 8 2 - 3(6 - tfe 2 ^ 2 ^ + 3(6' - c') 2 e 2 ( b - 2c - Q )l . (A10) 



Finally, the twice contracted Bianchi identities V^G* j = are 

Goo + (d + b + 2c)G 00 +aGu + b G 22 + c(G 33 + G 44 ) = e d - a [G' 01 + (b' + 2c! + 2d')G 01 ] , (All) 
G01 + (2d + b + 2c)G i = e d - a [G' n + d! G 00 + {b' + 2c' + d')G n - b' G 22 - c'(G 33 + G 44 )] . (A12) 



APPENDIX B: TANGHERLINI SOLUTION AND INVARIANTS 

Tangherlini's solution [l4[ in five dimensions reads 

ds 2 = (l dt 2 + (l - 5) ^ dr * + ^ (Bl) 

with G related to the ADM mass by 

3-7T 

M ADM = — C. (B2) 
o 

For the curvature invariants, we find n = 3C/r 4 , Pq = /3 2 = 70 = 72 = 6% = <5 2 = G/r 4 , /?i = 7 i = 0, and so the 
curvature invariants 7, I and I defined in appendix [A] are 

ri2 (~<2 (~i2 

7 = 36--, ^=44—, 1 = 72—. 

At the event horizon, the following combinations are particularly simple: Jn = Ir EH /36, K~n = Ir% H /4A, In = 
Ir^ H /72. Since r 2 EH = G, In = Jn = Kn = 1 when evaluated at the horizon. 
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